NASA Technical Memorandum 102733 


AEROELASTIC ANALYSIS OF WINGS USING THE EULER 
EQUATIONS WITH A DEFORMING MESH 


BRIAN A. ROBINSON 
JOHN T. BATINA 
HENRY T. Y. YANG 


(NASA-TM-102713) 
WINGS USING THt 

deforming MESH 


AEMULASTIC ANALYSIS CF 
EiULER bQUATIUNS WITH A 
(NASA) 11 p CSCLOIA 


03 / 0 ? 


N? 1-1 00 


Unc 1 <js 
0310830 


NOVEMBER 1990 


NASA 

National Aeronautics and 
Space Administration 

Langley Research Center 

Hampton. Virginia 23665 



AEROELASTIC ANALYSIS OF WINGS USING THE EULER EQUATIONS 
WITH A DEFORMING MESH 


Brian A. Robinson* 
McDonnell Aircraft Company 
St. Louis, Missouri 63166 

John T. Batina** 

NASA Langley Research Center 
Hampton, Virginia 23665-5225 

Henry T. Y. Yangt 
Purdue University 
West Lafayette, Indiana 47907 


Abstract 

Modifications to the CFL3D three-dimensional unsteady 
Euler/Navier-Stokes code for the aeroelastic analysis of 
wings are described. The modifications involve including a 
deforming mesh capability which can move the mesh to 
continuously conform to the instantaneous shape of the 
aeroelastically deforming wing, and including the structural 
equations of motion for their simultaneous time-integration 
with the governing flow equations. Calculations were 
performed using the Euler equations to verify the 
modifications to the code and as a first-step toward 
aeroelastic analysis using the Navier-Stokes equations. 
Results are presented for the NACA 0012 airfoil and a 45 
sweptback wing to demonstrate applications of CFL3D for 
generalized force computations and aeroelastic analysis. 
Comparisons are made with published Euler results for the 
NACA 0012 airfoil and with experimental flutter data for the 
45° sweptback wing to assess the accuracy of the present 
capability. These comparisons show good agreement and, 
thus, the CFL3D code may be used with confidence for 
aeroelastic analysis of wings. The paper describes the 
modifications that were made to the code and presents results 
and comparisons which assess the capability. 


Nomenclature 

nondimensional distance from midchord to elastic 
axis 


k| generalized stiffness of mode i 
m airfoil mass per unit span 
m j generalized mass of mode i 
freestream Mach number 
qj generalized displacement of mode i 

Q generalized force in mode i 

Q dynamic pressure, j u ~ 

5 nondimensional dynamic pressure, (U^/(ba) a VJI))^ 

r a airfoil radius of gyration about elastic axis 

t time, nondimensionalized by freestream speed of 

sound and airfoil chord, Ta^c 

T time, seconds 

Uj load vector 

IL streamwise freestream speed 


su freestream speed of sound 

A r generalized aerodynamic force resulting from 

^ pressure induced by mode j acting through mode i 

b semichord, c/2 

c reference chord 

Cj generalized damping of mode i 

Cp pressure coefficient 

k reduced frequency, <oc/2U.. 
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X| state vector 

x a nondimensional distance from elastic axis 
to mass center 

mean angle of attack 

e Integral of state-transition matrix 

p, mass ratio 

p_ freestream density 

a real part of Laplace transform variable 
<t> state-transition matrix 

to angular frequency 

o) h uncoupled natural frequency of bending mode 

co uncoupled natural frequency of torsion mode 


Introduction 

In recent years, there has been increased interest in the 
development of aeroelastic analysis methods involving 
computational fluid dynamics techniques. 1 This research has 
been highfy focused on developing finite-difference codes for 
the solution of the transonic small-disturbance 2 * 3 and full 
potential 4 * 5 equations, although efforts are currently 
underway at the higher equation levels as well. 6 ’ 11 For 
example, Bendlksen and Kousen 6 presented transonic flutter 
results for two-degree-of-freedom (plunging and pitching) 
airfoils by simultaneously integrating the structural 
equations of motion with the two-dimensional unsteady Euler 
equations. The Euler equations were Integrated using a 
Runge-Kutta time-stepping scheme involving a finite- 
volume spatial discretization and a moving mesh. The 
instantaneous mesh was taken to be a linear combination of 
meshes corresponding to rigid plunging and pitching of the 
airfoil, in a following study, Kousen and Bendiksen 7 applied 
their method of Ref. 6 to investigate the nonlinear 
aeroelastic behavior of two-degree-of-freedom airfoils at 
transonic speeds. In Ref. 7, transonic flutter instabilities 
were shown to lead to stable limit-cycle oscillations. Wu, 
Kaza, and Sankar 8 have time-integrated the unsteady 
compressible Navier-Stokes equations for airfoils 
undergoing one- and two-degree-of-freedom aeroelastic 
motions. In Ref. 8, flutter characteristics of airfoils at high 
angles of attack were investigated including cases with stall 
flutter. The method of Ref. 8 has also been recently applied 
by Reddy, Srivastava, and Kaza 9 to study the effects of 
rotational flow, viscosity, thickness, and shape on the 
transonic flutter dip phenomena. The study concluded that 
the influence of these effects on flutter, for the cases 
considered, was small near the minimum of the flutter dip, 
but may be large away from the dip. Guruswamy 10 has 
demonstrated the simultaneous time integration of the three- 
dimensional Euler equations along with the structural 
equations of motion. The capability was demonstrated in a 
time-marching flutter analysis performed for a rectangular 
wing with a parabolic-arc airfoil section. Finally, Rausch, 
Batina, and Yang 1 1 have presented Euler aeroelastic results 
for two degree-of-freedom airfoils using a flow solver based 
on unstructured grids. A novel aspect of the capability of 
Ref. 1 1 is the dynamic mesh algorithm which is used to move 
the mesh so it continuously conforms to the instantaneous 
position of the airfoil. The algorithm is completely general 
in that it can treat realistic motions and deformations of 
multiple two-dimensional bodies. 

Many of the methods that are currently being developed 
for aeroelastic analysis assume that the mesh moves rigidly 
or that the mesh shears as the body deforms. These 
assumptions consequently limit the applicability of the 
procedures to rigid-body or small-amplitude motions. 
These problems, for example, are easily demonstrated by 
considering an airfoil section of a wing whose aeroelastic 
deformation involves significant chordwise bending. As the 
airfoil section bends, grid lines that eminate from the 
concave surface of the section may collapse onto the 
neighboring grid lines. Similar difficulties can occur if the 
wing has significant spanwise bending. Also, for methods 
where the mesh is constrained such that airfoil sections of 
the wing can only pitch and plunge, chordwise deformation 
cannot be modeled and wing tip deflections must remain 
small. Therefore the purpose of the paper is to describe the 
implementation of a deforming mesh capability which 
effectively removes the rigid-body and small-amplitude 
limitations of previous methods. This capability was 
developed within the CFL3D unsteady Euler/Navier-Stokes 
code. 12 - 13 The deforming mesh capability is a general 
procedure, based on the dynamic mesh algorithm of Refs. 1 1 
and 14, which can move the mesh for realistic motions and 
structural deformations of wings. In addition, the structural 
equations of motion have been implemented within CFL3D to 


allow simultaneous time-integration with the governing flow 
equations for aeroelastic analysis. Calculations were 
performed using the Euler equations to verify the 
modifications to the code and as a first step toward 
aeroelastic analysis using the Navier-Stokes equations. 
Results are presented for the NACA 0012 airfoil and a 45° 
sweptback wing to demonstrate applications of CFL3D for 
generalized force computation and aeroelastic analysis. 
Comparisons are made with the Euler results of Ref. 1 1 for 
the NACA 0012 airfoil and with the experimental flutter data 
of Ref. 1 5 for the 45° sweptback wing to assess the accuracy 
of the present procedures. The paper describes the 
modifications that were made to the CFL3D code and presents 
results and comparisons which assess the capability. 


Euler Solution Algorithm 

In the present study the flow was assumed to be governed 
by the time-dependent Euler equations which may be written 
in conservation form as 


+ 

at ^ an ac 


=o 


(i) 


where the vector 6 represents the conserved variables 

divided by the Jacobian and £ , 6, A are the inviscid fluxes 
which have been transformed from the cartesian (x, y, z) 

coordinate system to generalized (£, n>0 coordinates. 
Equations (1) are solved within the CFL3D code by a three- 
factor, Implicit, finite-volume algorithm based on upwind- 
biased spatial differencing. The upwind-biased differencing 
involves either flux-vector splitting or flux-difference 
splitting implemented as a cell-centered discretization. 
Flux-limiting may also be used in the spatial differencing to 
determine values of the flow variables on the cell faces. For 
unsteady applications, the algorithm includes the grid speed 
metric terms that are necessary for time-accuracy with 
moving meshes, although the original scheme was limited to 
cases involving rigid-body plunge or pitch where the mesh 
moves without deformation. Modifications to the algorithm 
to include the terms arising from a deforming mesh, which 
are required for aeroelastic analysis, are described in the 
following section. 


Deforming Mesh Algorithm 

A deforming mesh algorithm was developed and 
implemented to move the mesh so that it continuously 
conforms to the instantaneous shape of the aeroelasticaily 
deforming wing. The method, based on that of Ref. 14, 
models the mesh as a spring network where each edge of each 
hexahedral cell represents a linear spring. The stiffness of 
each spring is inversely proportional to a specified power of 
the length of the edge. For example, along an edge 

( i ) - (i +1) , the stiffness k is 


k , + 1 “ 1 - 0 / [( x ,.,- x l ) 2 + (y t+ ,-y,) 2 + ( z i.,- 2 .) 2 ] P/ * 

2 

(2) 

In addition to the cell edges, springs are also placed along the 
diagonals of each cell face to control cell shearing, with 
spring stiffness defined similar to Eq. (2). The power p 
which appears in Eq. (2) is used to control the stiffness of 
the cells near the wing. These cells are typically very small 
in comparison with cells in the far field, and as such it has 
been found advantageous to increase the stiffness of the near 
field ceils to avoid excessive mesh distortion in this region. 
The stiffness is increased by increasing the value of p. In 
this study, p was typically set equal to two or three. 
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The procedure to move the mesh is described as follows. 
At each time step, the instantaneous positions of the points on 
the wing are prescribed while the points on the outer 
boundary are held fixed. The displacements of the interior 

points 5 x ,5 y ,5 z are determined by solving the static 

equilibrium equations which result from a summation of 
forces at each point in the x, y, and z directions. This 
solution is approximated by using a predictor-corrector 
procedure which first predicts the displacements of the 
interior points by a linear extrapolation of displacements 
from the two previous time levels according to 


= 2 5 w 


- 6 . 


n-1 
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(a) original grid (for reference). 
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and then corrects the displacements using several Jacobi 
iterations of the static equations written as 


(4a) 
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For the applications performed in the present study, two to 
four Jacobi iterations were sufficient to accurately move the 
mesh. 

To demonstrate mesh movement using the deforming 
mesh algorithm, consider the coarse grid about a NACA 0012 
airfoil that is shown in Fig. 1 (a). The grid is of C-type 
mesh topology and has 43 points in the •wrap-around* 
direction and 11 points in the outward direction. It is used 
only to illustrate how the mesh moves. In this example, the 
airfoil was plunged for one cycle of sinusoidal motion with an 
amplitude of one chordlength. The mesh at the maximum 
plunge displacement is shown in Fig. 1 (b) and the mesh at 
the minimum plunge displacement is shown in Fig. 1 (c). 
The mesh moves smoothly as the airfoil plunges, and the 
procedure is completefy general in that it can treat realistic 
airfoil or wing motions including aeroelastic transient- type 
motion. 

Since the mesh can now deform to accommodate the 
aeroelastically deforming wing, the flow solver in the CFL3D 
code was modified to include an extra term in the time- 
discretization of the governing equations, to account for the 
mesh deformation. Specifically, the modification involves 
the change in cell volume when the mesh deforms. The 

algorithm changes are derived by first writing a£y at in the 




(c) plunge downward one chordlength. 

Fig. 1 Sequence of grids about the NACA 0012 airfoil 
which illustrates how the mesh moves for a 
plunging airfoil. 
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governing equations (Eq. (1)) as 3(QV)/ at, where V is the 
cell volume. The original flow solver assumed that the cell 
volume does not change in time so that 

a(Q v) _ v ao 

at at (5a) 

However, if the cell volume changes in time, as it does when 
the mesh deforms, the time derivative becomes 


WY1- V i2 + Q av 

at at at 


(5b) 


which requires the implementation of the Qav/at term 
within the CFL3D algorithm. 


Pulse Transfer-Functio n Analysis 
Generalized aerodynamic forces that are used in 
aeroelastic analyses are typically obtained by calculating 
several cycles of a harmonically forced oscillation with the 
determination of the forces based on the last cycle of motion. 
This method of harmonic oscillation requires one flow 
calculation for each value of reduced frequency that is of 
interest. In contrast, the generalized forces may be 
determined for a wide range of reduced frequency in a single 
flow calculation by the pulse transfer-function analysis. In 
the pulse analysis, the forces are computed indirectly from 
the response of the flow field due to an exponentially shaped 
pulse. The analysis assumes that the system is linear which 
is a reasonable assumption even for transonic cases, since 
experience has shown that the response for harmonic or 
aeroelastic motion is, in general, locally linear for small 
amplitudes of oscillation. 


Time-Marchinfl Aeroelastic Analysis 
In this section the aeroelastic equations of motion, the 
time-marching solution procedure, and the modal 
identification technique are described. 

Aeroelastic Equations of Motion 

The aeroelastic equations of motion that were 
incorporated within CFL3D were derived by assuming that 
the general motion of the wing is described by a separation of 
time and space variables in a finite modal series. 16 This 
modal series involves the summation of free vibration modes 
weighted by generalized displacements. Considering 
Lagrange's equations leads to the equations of motion which 
can be written for each mode i as 

ma+cd+ka-Q 

1 M i H l i (6) 

where qj is the generalized normal mode displacement, mj is 
the generalized mass, Cj is the generalized damping, k; is the 
generalized stiffness, and Qj is the generalized force 
computed by integrating the pressure weighted by the mode 
shapes. 

Time-Marchina Solution 

The aeroelastic equations of motion are integrated in 
time in a manner similar to that described by Edwards, et 
al .17,18 The formulation is implemented herein for 
multiple degrees-of-freedom or mode shapes of a wing 
following Ref. 16. Each normal mode equation represented 
by Eq. (6) can be expressed in state-space form as 


Xj = A x i + Bu . 


(7) 


where A and B are coefficient matrices that result from the 
change of variables Xj - [qj qj] T and Uj is the nondimensional 
generalized force weighted by mode i. Equation (7) is 
integrated in time using the modified state-transition matrix 
structural integrator 18 implemented as a predictor- 
corrector procedure, which first uses a linear extrapolation 
of Uj from previous time steps as 

8, +1 =‘*>x%eB(3u^-u^ 1 )/2 (8 a) 

n+1 n +1 „rv+1 

to compute * } , the prediction for x t . Then, * is 

used to compute the flow field and evaluate the load vector 
n+1 

°j . These values are then used in the corrector step to 
determine xj n+1 given by 


(8b) 

In Eqs. (8a) and (8b), <t> is the state-transition matrix and 
9 is the integral of the state transition matrix from time 
step n to n+1 . 

Modal Identification Technique 

Damping and frequency characteristics of the aeroelastic 
responses are estimated from the response curves by using 
the modal identification technique of Bennett and 
Desmarais. 19 The modal estimates are determined by a least 
squares curve fit of the responses of the form 



i - 1,2. ... 

where m is the number of modes. 


Results and Discussion 

Results are presented in this section for the NAGA 0012 
airfoil and a 45° sweplback wing, computed using CFL3D, to 
verify the deforming mesh capability and to assess the code 
for aeroelastic analysis. The accuracy of these results Is 
determined by making detailed comparisons with CFL3D 
calculations performed using a rigidly moving mesh, 
published results obtained using alternative computational 
methods, and available experimental data. 



Fig. 2 Partial view of 159 x 49 C-type mesh about the 
NACA 0012 airfoil. 
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NACA 0012 Airfoil Results 

Calculations were performed for the NACA 0012 airfoil 
by using the CFL3D code run in a 2-D mode. The results 
were obtained using a 159 x 49 C-type mesh, a partial view 
of which is shown in Fig. 2. The outer boundaries of the 
mesh were located approximately fifteen chordlengths from 
the airfoil and there are 110 points which lie on the airfoil 
surface. Calculations were performed for the airfoil at - 
0.8 and zero degrees angle of attack. In these calculations, 
the Euler equations were solved using the CFL3D code with a 
third-order accurate upwind-biased spatial discretization 
and flux-vector splitting. Steady results are compared with 
the experimental pressure data of Ref. 20 and the unsteady 
results are compared with the parallel computational results 
of Ref. 11, obtained using the 2-D Euler code of Ref. 21. 

Steady Pressure Comparisons^ - The calculated steady 
pressure distribution along the upper surface of the airfoil 
is compared with the experimental data in Fig. 3. These 
pressures indicate that there is a moderately strong shock 
wave near the airfoil midchord, which is accurately 
predicted by CFL3D in both strength and location. Such a 
good comparison suggests that viscous effects for this case 
are relatively small and, thus, the flow can be modeled 
accurately by the Euler equations. 

npnfiralized Force C omparisons^ - Generalized 
aerodynamic forces for the NACA 0012 airfoil are presented 
in Fig. 4. The results are plotted as real and imaginary 
components of the unsteady forces, A,j f as a function of 
reduced frequency k. Both plunge and pitch-about-the- 
quarter-chord motions were considered, which are defined 
as modes 1 and 2, respectively. Thus, for example, A^ 
the lift coefficient due to pitching. Figure 4 compares 
results obtained using several different methods including: 
(1) the pulse analysis with a rigidly moving mesh; (2) the 
pulse analysis with a deforming mesh; (3) harmonic motion 
with a deforming mesh; and (4) the harmonic Euler results 
from Ref. 11. With the rigidly moving mesh, the mesh 
simply translates for airfoil plunge motion and rotates for 
airfoil pitch motion. The harmonic results were obtained at 
six values of reduced frequency: k - 0.0, 0.125, 0.25, 0.5, 
0.75, and 1.0. The amplitudes of motion were 0.01 
chordlengths and 0.1 degrees for plunge and pitch, 
respectively, in both harmonic and pulse analyses. 

As shown in Fig. 4, the forces from the pulse analysis 
obtained using the deforming mesh agree very well with the 
forces obtained using the rigidly moving mesh. This good 
agreement between the two sets of forces tends to verify the 
deforming mesh capability that was implemented within 
CFL3D. As further shown in Fig. 4, the pulse results agree 
well with the forces from the harmonic analysis, for both 
plunge and pitch motions, for the entire range of reduced 
frequency that was considered. The harmonic analysis, 
however, is considered to be the more accurate of the two 
sets of calculations, since the local linearity assumption in 
the pulse analysis is questionable for transonic cases. 
Furthermore, the generalized forces determined using the 
harmonic analysis agree well with the Euler forces of 
Rausch, et al. 11 which gives additional confidence in the 
accuracy of the deforming mesh capability that was 
implemented. 

ARroelastic Comparisons. - Aeroelastic results are 
presented for a much-studied case designated as Case A of 
Isogai, 22 which has normal modes similar to those of a 
stream wise section near the tip of a sweptback wing. The 
wind-off bending and torsion natural frequencies are 71.33 
and 535.65 rad/sec, respectively. The pivot point for the 
bending mode is located 1 .44 chordlengths upstream of the 
leading edge of the airfoil. The pivot point for the torsion 
mode is 0.068 chordlengths forward of midchord. These 



Fig. 3 Comparison of upper surface steady pressure 
distributions on the NACA 0012 airfoil at - 0.8 

and c X q «* 0°. 

Pulse - rigid mesh 

Pulse - deforming mesh 




Fig. 4 Comparisons of generalized aerodynamic forces for 
the NACA 0012 airfoil at - 0.8 and a 0 - 0°. 
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Table 1 Comparisons between aeroelastic solutions for the NACA 0012 airfoil 
at - 0.8 and ao - 0° for Case A. 


<3 

Method 

Model 

Mode 2 1 

a/(0 a 

0)/<i)a 

C/C0a 

<o/(Da 

.2 

CFL3D time-marching 
Euler (Rausch) time-marching 

-.011 

-.011 

.794 

.790 

-.091 

-.068 

5.363 

5.353 

.5 

CFL3D time-marching 
Euler (Rausch) time-marching 

.004 

.000 

.914 

.913 

-.185 

-.148 

5.347 

5.349 

.8 

CFL3D time-marching 
Euler (Rausch) time-marching 

.026 

.017 

1.027 

1.022 

-.173 

-.223 

5.270 

5.317 


X 10* 3 
1 i- 


CFL3D 

Euler (Rausch) 


q 2 o 


-i *- 
1 1- 


.; lvAA : 


0 - 0.2 

5-0.5 



Fig. 5 Comparisons of generalized displacements for the 
NACA 0012 airfoil at - 0.8 and a 0 - 0° for 


Case A. 


mode shapes and natural frequencies were determined by 
performing a free vibration analysis with the aeroelastic 
equations written in the traditional form of plunge and pitch 
degrees of freedom. In this analysis, the following 
structural parameter values were used: a - -2.0, x a - 1.8, 
r a « 1.865, <o h « 100 rad/sec, and co a - 100 rad/sec. Also, 
the airfoil mass ratio was p - 60. Generalized displacements 
corresponding to the bending and torsion modes are defined as 
q^ and q 2 , respectively. Initial conditions for the time- 
marching aeroelastic analysis were i (0) - 2.0 and q 2 ( 0 ) * 
0 . 01 . 

Aeroelastic results for Case A were obtained for several 
values of nondimensional dynamic pressure including Q » 
0.2, 0.3, 0.4, 0.5, 0.6, 0.7, and 0.8, to obtain conditions 
which bracket the flutter point. Figure 5 shows time 
responses of generalized displacement of the second coupled 
mode for 3 - 0.2, 0.5, and 0.8 which correspond to stable, 
near neutrally stable, and unstable aeroelastic behavior, 
respectively. Also plotted are the corresponding responses 
reported in Ref. 11. A comparison of these responses 
indicates that the time-marching aeroelastic results from 
the CFL3D code agree well with those from the Euler code of 
Ref. 11, which tends to verify the aeroelastic modeling 
procedures that were implemented. Shown In Fig. 6 are the 
two-mode curve fits of the CFL3D responses which are 
excellent approximations to the original data. The component 
modes from these curve fits are shown in Fig. 7 for the three 


X 10- 3 
1 i- 


Data 

Two-mode fit 


q 2 




0 « 0.2 



Amplitude 


Mode 1 


Mode 2 

x 10- 3 
1 r 


-1 


G - 0.2 


r Amplitude 


Amplitude 


0 

-tLV w v a - o.5 



Fig. 6 Effects of nondimensional dynamic pressure on the 
generalized displacement of the second coupled 
mode for the NACA 0012 airfoil at M„ » 0.8 and 

olq - 0° for Case A. 


Fig. 7 Effects of nondimensional dynamic pressure on the 
component modes of the aeroelastic system for the 
NACA 0012 airfoil at M_ * 0.8 and a 0 - 0° for 

Case A. 
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values of 0 that were considered in Figs. 5 and 6. The results 
of Fig. 7 show that the component modes consist of a 
dominant mode corresponding to bending (mode 1) and a 
second higher-frequency mode corresponding to torsion 
(mode 2). Damping and frequency estimates from this 
analysis are compared with similar values from Rausch, et 
ai.11 in Table 1. These comparisons Indicate that the CFL3D 
values correlate well with the Euler values from Ref. 1 1 . 
Also, the flutter value for G computed by quadratic 
interpolation of the damping values, was 0.48 for CFL3D 
which compares with 0.50 as reported in Ref. 11. Linear 
theory at M« - 0.8, which of course does not include 
transonic effects, predicts a much higher flutter value of 
1.89. 


Calculations were performed for a simple well-defined 
wing, to assess the CFL3D code for three-dimensional 
aeroelastic applications. The wing that was analyzed was a 
semispan wind-tunnel-wall-mounted model that has a 
quarter-chord sweep angle of 45°, a panel aspect ratio of 
1.65, and a taper ratio of 0.66. 15 The wing is an AGARD 
standard aeroelastic configuration which was tested in the 
Transonic Dynamics Tunnel (TDT) at NASA Langley Research 
Center. A planview of the wing is shown in Fig. 8. The wing 
has a NACA 65A004 airfoil section and was constructed of 
laminated mahogany. In order to obtain flutter for a wide 
range of Mach number and density conditions in the TDT, 
holes were drilled through the wing to reduce its stiffness. 
To maintain the aerodynamic shape of the wing, the holes 
were filled with a rigid foam plastic. A photograph of the 



MttfcM. -160Hz 


l * f 2 - 3117 Hz 



Mod* 3,(3 -48.35Hz 


* 4. 1 4 -9134 Hz 


Fig. 10 Oblique projections of natural vibration modes of 
45° sweptback wing. 




1 2J 2 -38.17 Hz 
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Fia. 8 Planview of 45° sweptback wing. 


Fig. 11 Deflection contours of natural vibration modes of 
45° sweptback wing. 




45° sweptback wing in the NASA Langley Transonic 
Dynamics Tunnel. 


Fig. 12 Partial view of 193 x 33 x 41 C-H-type mesh 
about the 45° sweptback wing. 


ORIGINAL PAGE 

black and v/hitb photograph 






wing mounted in the TDT is shown in Fig. 9. The wing is 
modeled structurally using the first four natural vibration 
modes which are illustrated in Figs. 10 and 11. Figure 10 
shows oblique projections of the natural modes while Fig. 11 
shows the corresponding deflection contours. These modes 
which are numbered 1 through 4 represent first bending, 
first torsion, second bending, and second torsion, 
respectively, as determined by a finite element analysis. 
The modes have natural frequencies which range from 9.6 Hz 
for the first bending mode to 91 .54 Hz for the second torsion 
mode. 

Aeroelastic results were obtained for the 45° sweptback 
wing using a 193 x 33 x 41 C-H-type mesh, a partial view 
of which Is shown in Fig. 12. Calculations were performed 

for the wing at -0.9 and zero degrees angle of attack. 
In these calculations, the Euler equations were solved using 
the CFL3D code with a second-order accurate upwind-biased 
spatial discretization and flux-vector splitting. Aeroelastic 
transients were obtained for several values of dynamic 
pressure Q, to obtain conditions which bracket the flutter 
point. Figure 13 shows time responses of generalized 
displacement of the first bending mode for Q-0.9 Q exp. 1 -0 
Q exp. arid 1.1 Q exp, where O exp is the experimental 
flutter dynamic pressure value. Also shown in Fig. 13 are 
the two-mode curve fits of the responses which are very 
good approximations to the original data. The component 
modes from these curve fits are shown in Fig. 14 for the 
three values of Q that were considered in Fig. 13. The 
results of Fig. 14 show that the component modes consist of a 
dominant mode corresponding to first bending (mode 1) and a 
higher-frequency damped mode corresponding to first 
torsion (mode 2). Also, the flutter value for O computed by 
quadratic interpolation of the damping values was 0.92 
Qexp- Although the calculated value is slightly low in 
comparison with the experimental value, it is within 1% of 
the value predicted by the CAP-TSD transonic small- 
disturbance code, 1 ® which tends to verify the computational 
aeroelasticity methods of the present study. 


Concluding Remarks 

Modifications to the CFL3D three-dimensional unsteady 
Euler/Navier-Stokes code for the aeroelastic analysis of 
wings were described. The modifications involve including a 
deforming mesh capability which can move the mesh to 
continuously conform to the instantaneous shape of the 
aeroelastically deforming wing, and including the structural 
equations of motion for their simultaneous time-integration 
with the governing flow equations. Calculations were 
performed using the Euler equations to verify the 
modifications to the code and as a first step toward 
aeroelastic analysis using the Navier-Stokes equations. 

Results were presented for the NACA 0012 airfoil and a 
45° sweptback wing to demonstrate applications of CFL3D 
for generalized force computations and aeroelastic analysis. 
Detailed comparisons were made with published Euler 
results for the NACA 0012 airfoil which indicated very good 
agreement for generalized forces due to harmonic motion in 
pitch or plunge, and good agreement for aeroelastic 
transients corresponding to stable, neutrally stable, and 
unstable aeroelastic behavior. This favorable correlation 
tends to verify the deforming mesh capability and the 
aeroelastic modeling procedures that were implemented 
within CFL3D. Aeroelastic transients were obtained for a 
45° sweptback wing which also demonstrated stable, 
neutrally stable, and unstable behavior. The resulting 
flutter dynamic pressure, determined by interpolation of the 
dominant damping values, was within 6% of the 
experimental flutter value and within 1% of a transonic 
small-disturbance result. 
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